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ABSTRACT 

We introduce project NIHAO (Numerical Investigation of a Hundred Astrophysical 
Objects), a set of 100 cosmological zoom-in hydrodynamical simulations performed 
using the GASOLINE code, with an improved implementation of the SPH algorithm. 
The haloes in our study range from dwarf (M 200 5 x IO^Mq) to Milky Way 

(M200 ^ 2 X IO^^Mq) masses, and represent an unbiased sampling of merger his¬ 
tories, concentrations and spin parameters. The particle masses and force softenings 
are chosen to resolve the mass profile to below 1% of the virial radius at all masses, 
ensuring that galaxy half-light radii are well resolved. Using the same treatment of 
star formation and stellar feedback for every object, the simulated galaxies reproduce 
the observed inefficiency of galaxy formation across cosmic time as expressed through 
the stellar mass vs halo mass relation, and the star formation rate vs stellar mass 
relation. We thus conclude that stellar feedback is the chief piece of physics required 
to limit the efficiency of star formation in galaxies less massive than the Milky Way. 

Key words: galaxies: evolution - galaxies: formation - galaxies: dwarf - galaxies: 
spiral - methods: numerical - cosmology: theory 


1 INTRODUCTION 

Understanding the formation of spiral galaxies (including 
our own Galaxy - the Milky Way) has been at the fore¬ 
front of theoretical astrophysics for decades. While there is 
a broadly-sketched paradigm for the formation of disk galax¬ 
ies (White & Rees 1978; Fall & Efstathiou 1980; Mo, Mao, 
& White 1998; Dutton et al. 2007), the details have yet 
to be understood. The simulation of realistic spiral galaxies 
from cosmological initial conditions has been a formidable 
challenge. Early (low-resolution) simulations were plagued 
by over-cooling, and angular momentum losses resulting in 
compact bulge dominated galaxies (e.g., Navarro & Stein- 
metz 2000). Recently, it has been shown that higher spatial 
resolution coupled to more realistic models for star forma¬ 
tion and stellar feedback can solve these problems, result¬ 
ing in spiral galaxies with realistic sizes and bulge fractions 
(Guedes et al. 2011; Brook et al. 2012; Aumer et al. 2013; 
Marinacci et al. 2014; Grain et al. 2015). 

* lwang@mpia.de 


One of the most fundamental constraints for galaxy for¬ 
mation models is the relation between the stellar mass of a 
galaxy, Mstar, and the mass of its host dark matter halo, 
Afhaio- This relation encodes the global efficiency at which 
galaxies turn gas into stars. In the context of AGDM, a 
model that matches the Mstar — Mhaio relation will also re¬ 
produce the observed galaxy stellar mass function, which 
has been a benchmark for galaxy formation models for over 
a decade (e.g., White & Frenk 1991, Benson et al. 2003; 
Crain et al. 2009; Oppenheimer et al. 2010). An advantage 
of the Mstar vs Mhaio relation is that single galaxies can be 
compared, whereas comparing with the stellar mass function 
requires a cosmologically representative volume. 

In recent years the Mstar vs Mhaio relation has been 
studied extensively using observations of weak gravita¬ 
tional leasing, satellite kinematics, and halo abundance 
matching (e.g., Yang et al. 2003, 2012; Mandelbaum 
et al. 2006; Gonroy & Wechsler 2009; Moster et al. 2010; 
More et al. 2011; Leauthaud et al. 2012; Behroozi et al. 2013; 
Garrison-Kimmel et al. 2014; Kravtsov et al. 2014; Hud¬ 
son et al. 2015). One of the key results of these studies is 
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Figure 1. Concentration vs halo mass (left) and spin parameter vs halo mass (right) for the sample of haloes presented in this paper. 
Concentrations are from NFW fits. Solid and dashed lines show the mean and Itr scatter from the low-resolution parent simulations in 
the Planck cosmology from Dutton & Maccio (2014). 


that star formation is inefficient - the maximum efficiency 
is ~ 25%, and independent of redshift, occurring at a halo 
mass similar to that of the Milky Way (i.e., ~ lO^^M©). At 
both higher and lower halo masses the efficiency drops even 
further. 

Reproducing this relation both at redshift z = 0 and 
earlier times, in a fully cosmological context has been a 
challenge, with many simulations and semi-analytic models 
predicting either too many stars at redshift 2 = 0 , or ear¬ 
lier times (Weinmann et al. 2012). Even some of the largest 
zoom-in and cosmological volume simulations do not repro¬ 
duce the observed 2 = 0 relation between stellar and halo 
mass (e.g., Guedes et al. 2011; Marinnaci et al. 2014; Vo- 
gelsberger et al. 2014). Only a handful of cosmological sim¬ 
ulations have been able to match the stellar vs halo mass 
relation over a wide range of halo masses from dwarfs to the 
Milky Way (Di Cintio et al. 2014; Hopkins et al. 2014; Schaye 
et al. 2015). On the other hand these simulations were sub¬ 
ject to either limited resolution on disc galaxy scales (Schaye 
et al. 2015) or based on a handful of objects with unclear 
selection function (DiCintio et al. 2014; Hopkins et al. 2014) 
and thus not able to recover the variety of formation histo¬ 
ries of real galaxies. 

The next step is then to create a large sample of galaxies 
with very high resolution over a range of halo masses with 
an unbiased sampling of the mass accretion histories. The 
NIHAC0 (Numerical Investigation of a Hundred Astrophys- 
ical Objects) project aims to achieve these goals with ~ 100 
hydrodynamical cosmological zoom-in simulations with halo 
masses ranging from M^aio ~ 5 x 10® to ~ 2 x lO^^M©. 

In this paper we describe the numerical methods and 
sample selection and present results on the growth of stellar 
mass across cosmic time from the first 88 NIHAO simula- 

^ Nihao is the Chinese word for “hello”. 


tions. Subsequent papers will discuss: dark halo shapes (But- 
sky et al. 2015); central dark matter density slopes (Toilet 
et al. 2015); Tully-Fisher relations (Dutton et al. in prep.); 
and gas content (Stinson et al. 2015). This paper is organized 
as follows: The simulations including sample selection, hy¬ 
drodynamics, star formation and feedback are described in 
^ results including the stellar mass vs halo mass and star 
formation rate vs halo mass relation are presented in ^and 
a summary in m 


2 SIMULATIONS 

The simulations presented here are a series of fully cosmolog¬ 
ical “zoom-in” simulations of galaxy formation run in a flat 
ACDM cosmology with parameters from the Planck Collab¬ 
oration et al. (2014): Hubble parameter Ho= 67.1 kms“^ 
Mpc“^, matter density flm = 0.3175, dark energy density 
Ha = 1 — Hm — Hr = 0.6824, radiation density Hr = 0.00008, 
baryon density Hb = 0.0490, power spectrum normalization 
as = 0.8344, power spectrum slope n = 0.9624. 

2.1 Sample selection and Initial Conditions 

The simulated galaxies were selected from two cosmological 
boxes of sides 60 and 20 Mpc//i from Dutton & Maccio 
(2014), as well as a new simulation of size 15 Mpc//i 
with 400® particles. We select halo masses from 9.5 < 
log^Q(M 2 oo/M 0 ) < 12.3. Haloes are selected without con¬ 
sidering the halo structure or merger history. For technical 
reasons with the zoom-in technique we require the haloes to 
be “isolated”. Formally, we only consider haloes that have 
no other haloes with mass greater than one-fifth of the virial 
mass within 3 virial radii. 95% of all parent haloes satisfy 
this constraint, so it does not bias our sample in any sig¬ 
nificant way. The concentration vs mass and spin vs mass 
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Table 1. Dark and gas particle masses, m, and force softenings, e, for our zoom-in simulations. 


box.level halo mass range 
(Mq) 

nigas 

(Me) 

^gas 

(pc) 

muM 

(M©) 

Com 

(pc) 

60.1 

6.96 X 10^1 - 2.79 x lO^^ 

3.166x10® 

397.9 

1.735x10® 

931.4 

60.2 

8.89 X lO^o - 5.55 x 10“ 

3.958 xlO** 

199.0 

2.169x10® 

465.7 

60.3 

2.63 X 10i° - 7.12 X lOi® 

1.173x10"^ 

132.6 

6.426x10"^ 

310.5 

20.1 

4.99 X 10® - 2.94 x lO^® 

3.475x10® 

88.4 

1.904x10'^ 

207.0 

15.1 

5.41 X 10® - 9.26 X 10® 

2.087x10® 

74.6 

1.144x10"^ 

174.6 

20.2 

4.36 X 10® 

1.466x10® 

66.4 

8.032x10® 

155.2 

15.2 

3.54 X 10® 

6.184x10® 

49.7 

3.389x10® 

116.4 


relations from dark matter only simulations are shown in 
Fig.m The solid and dashed lines show the mean and scatter 
from the parent simulations from Dutton & Maccio (2014) 
showing that our selected haloes span a wide range in con¬ 
centrations and spins. 

Initial conditions for zoom-in simulations were cre¬ 
ated using a modified version of the GRAFIC2 package 
(Bertschinger 2001) as described in Penzo et al. (2014). We 
chose the refinement level in order to maintain a roughly 
constant relative resolution (i.e., £DM/I?vir ~ 0.003), and 
~ 10® dark matter particles per halo. This allows us to reli¬ 
ably probe the mass profile down to 1% of the virial radius 
or better across the full range of halo masses we simulate. 
The motivation for achieving a constant relative resolution, 
rather than a fixed physical resolution is that we wish to 
resolve the dark matter mass prohle down to the galaxy 
half-light radius, which is typically 1-2% of the virial radius 
(Kravtsov 2013). 

Fig. [2] shows the mass resolution of our simulations 
compared to a number of state-of-the-art zoom-in (ERIS - 
Guedes et al. 2011; Governato et al. 2012; Aumer et al. 2013; 
Marinacci et al. 2014; FIRE - Hopkins et al. 2014; Sawala 
et al. 2015) and large volume (ILLUSTRIS - Vogelsberger 
et al. 2014; EAGLE - Schaye et al. 2015) simulations. NI- 
HAO is by far the largest set of zoom-in simulations resolved 
with ~ 10® particles per halo. The large volume simulations 
(EAGLE, ILLUSTRIS) contain many thousands of galaxies, 
but the fixed force and mass resolution means that only the 
highest mass haloes are well resolved. 

The force softenings and particle masses for the highest 
rehnement level for each simulation are given in Table [T] 
Note that the ratio between dark and gas particle masses 
is initially the same as the cosmological dark/baryon mass 
ratio DoM/Db ~ 5.48, and the gas (and star) particle force 
softenings are set to be (flDivi/Db) ~ 2.34 times smaller 
than the dark matter particle softenings. Each hydrody¬ 
namic simulation has a corresponding dark matter-only sim¬ 
ulation run at the same mass and force resolution to enable 
a study of the effects of galaxy formation on dark halo struc¬ 
ture (Dutton et al. in prep). 

2.2 Hydrodynamics 

We use a new version of the N-body SPH solver GASOLINE 
(Wadsley et al. 2004). A complete description of the new 
gasoline code including tests is given in Wadsley et al. (in 
prep.). In their paper describing superbubble feedback, 
Keller et al. (2014) described the updates they made to 



Dark Halo Mass 

Figure 2. Number of dark matter particles per halo vs halo mass 
for NIHAO simulations (black circles) together with state-of-the- 
art cosmological simulations in the literature. 


GASOLINE. While we do not use their superbubble feedback, 
we employ their modified version of hydrodynamics that re¬ 
duces the formation of blobs and improves mixing. We thus 
refer to our version of GASOLINE as ESP-GASOLINe2. The 
biggest differences to the hydrodynamics in ESF-GASOLINe2 
come from the small change Ritchie & Thomas (2001) pro¬ 
posed for calculating Pjp^. Ritchie & Thomas (2001) also 
proposed modifying the density calculation to use equal 
pressures, but we do not use those densities in the simu¬ 
lations described here. 

Diffusion of quantities like metals and thermal energy 
between particles has been implemented as described in 
Wadsley et al. (2008). Metal diffusion is used, but ther¬ 
mal diffusion is not used because it is incompatible with the 
blastwave feedback that delays cooling. Gasoline2 includes 
several other changes to the hydrodynamic calculation. The 
Saitoh & Makino (2009) timestep limiter was implemented 
so that cool particles behave correctly when a hot blastwave 
hits them. To avoid pair instabilities, ESF-GASOLINe2 uses 
the Wendland C2 function for its smoothing kernel (Dehnen 
& Aly 2012). The treatment of artificial viscosity has been 
modified to use the signal velocity as described in Price 
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Figure 3. Edge-on views of a subset of NIHAO galaxies after processing through the Monte Carlo radiative transfer code SUNRISE. 
Images are 50 kpc on a side. 


(2008). We also increase the number of neighbor particles 
used in the calculation of the smoothed hydrodynamic prop¬ 
erties from 32 to 50. 

Cooling via hydrogen, helium, and various metal-lines 
in a uniform ultraviolet ionizing background is included 
as described in Shen et al. (2010) and was calculated us¬ 
ing CLOUDY (version 07.02; Ferland et al. 1998) These cal¬ 
culations include photo ionization and heating from the 
Haardt & Madau (private communication) UV background 
and Compton cooling and range from 10 to 10^ K. In the 
dense interstellar medium gas the extragalactic UV field is 
a reasonable approximation, and thus we do not impose any 
shielding. 

2.3 Star Formation and Feedback 

Within the hydrodynamic simulations, gas is eligible to form 
stars according to the Kennicutt-Schmidt Law when it sat¬ 
isfies a temperature and density threshold. Our fiducial runs 
adopt T < 15000 K and nth > 10.3 cm“®. The stars feed 


energy back into the interstellar medium (ISM) gas through 
blast-wave supernova feedback (Stinson et al. 2006) and ion¬ 
izing feedback from massive stars prior to their explosion as 
supernovae, referred to as “early stellar feedback” (Stinson 
et al. 2013). 

In GASOLINE, as in Stinson et al. (2013), the pre-SN 
feedback consists of cesf =10% of the total stellar flux, 
2 X 10®'^ erg of thermal energy per Mq of the entire stel¬ 
lar population, being ejected from stars into surrounding 
gas. Because of the increased mixing in ESF-gasoline2, 
the simulations required more stellar feedback to have their 
star formation limited to the abundance matching value at 
the Milky Way scale. Thus, we set eESF = 13%, which gives a 
better match to Behroozi et al. (2013) abundance matching 
results. Radiative cooling is left on for the pre-SN feedback. 

In the second, supernova, epoch of feedback, stars of 
mass 8 M© < Mstar < 40 M© eject both energy and 
metals into the interstellar medium gas surrounding the 
region where they formed. Supernova feedback is imple¬ 
mented using the blastwave formalism described in Stinson 
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Figure 4. Face-on views of a subset of NIHAO galaxies after processing through the Monte Carlo radiative transfer code SUNRISE. Images 
are 50 kpc on a side. 


et al. (2006). Since the gas receiving the energy is dense, it 
would quickly be radiated away due to its efficient cooling. 
For this reason, cooling is delayed for particles inside the 
blast region for ~ 30 Myr. 

2.4 Derived galaxy and halo parameters 

In this section we describe how the various galaxy and halo 
parameters we use in the rest of the paper are defined and 
measured. 

Haloes in our zoom-in simulations were identified using 
the MPI-l-OpenMP hybrid halo finder AHI0 (Knollmann & 
Knebe 2009; Gill et al. 2004). AHF locates local over-densities 
in an adaptively smoothed density field as prospective halo 
centers. The virial masses of the haloes are defined as the 
masses within a sphere containing A = 200 times the cosmic 
critical matter density. The virial mass, and size are denoted 

^ http://popia.ft.uam.es/AMIGA 


M 200 and i? 200 - The mass in stars, Mstai, is measured within 
a sphere of radius, rg^i = 0.2i?200- The star formation rate, 
SFR, is measured as the mass of stars formed inside rgai 
over the preceding 100 Myr. These parameters of the main 
halo/galaxy are given in Table lAll The simulation ID is 
named after the redshift z = 0 halo mass (Mvir/[MQ/h]) 
from the low resolution dark matter-only simulation. The 
actual halo masses in the hydrodynamic simulation may dif¬ 
fer due to: baryonic mass loss; differential evolution; and a 
different halo mass definition. 

Edge-on and face-on images for a selection of 30 galax¬ 
ies at 2 = 0 are shown in Figs. [3]&[31 Each image is 50 kpc 
on a side and was created using the Monte Carlo radiative 
transfer code sunrise (Jonsson 2006). The image brightness 
and contrast are scaled using arcsinh as described in Lup- 
ton et al. (2004). The image orientation was determined us¬ 
ing the angular momentum vector of the stars inside 5 kpc. 
We see a variety of galaxy morphologies, sizes, colors, and 
bulge-to-disc ratios. The galaxies are ordered according to 
their halo masses. The variation in morphology follows the 
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increase of halo mass. The lowest mass galaxies are small 
blue dots. Higher mass halos gradually extend into peanut 
shapes before filling ont entire discs. Once the galaxies be¬ 
come discs, their stellar popnlations start to include more 
old, red stars, and exhibit stronger absorption from dust 
lanes. These discs then start to show extended merger struc¬ 
tures. In some cases the disc has disappeared and all that 
is left is a smaller, red spheroid that represents the remnant 
of a merger. Even though all the images are taken from the 
z — 0 outputs, the mass sequence seems to replicate the 
evolutionary sequence of galaxies. 

In addition to the main target halo each simulation con¬ 
tains several smaller haloes within the high-res region, which 
we include in our analysis. We verify with convergence tests 
(see 33] below) that these poorer resolved haloes have similar 
global parameters to their higher resolution counterparts. 

3 RESULTS 

3.1 Stellar mass vs halo mass 

Fig. [S] shows the relation between stellar mass and virial 
mass for our simulated central galaxies compared to re¬ 
sults from three different halo abundance matching mea¬ 
surements; Behroozi et al. (2013), Moster et al. (2013), and 
Kravtsov et al. (2014). We show all three relations to rep¬ 
resent the lingering uncertainty in the abundance match¬ 
ing and halo occupation techniques. Most of the difference 
is apparent at low stellar masses, Mstar ^ 10® M©. That 
is the lower mass limit for the Sloan Digital Sky Survey 
(SDSS, Blanton et al. 2005). The smaller, deeper survey, 
GAMA (Driver et al. 2012), has pushed the stellar mass 
function to lower masses, but only the Kravtsov line uses 
data from the most up-to-date GAMA survey. The lines at 
low masses are dashed indicating that they are extrapola¬ 
tions of the higher mass results. The Kravtsov relation also 
differs from the other two at high masses because of how it 
takes into account intra-cluster light, namely by using the 
Bernardi et al. (2013) stellar mass function that extends sur¬ 
face brightness proHles to larger radii. Our simulations agree 
well with each of these relations, especially the most up-to- 
date relation from Kravtsov et al. (2014) at high masses. At 
low masses, the simulations most closely follow the extrap¬ 
olation of the Moster et al. (2013), which has the steepest 
slope, Mst ar oc M|qq. 

An additional systematic uncertainty to consider is due 
to form of the stellar initial mass function (IMF). The abun¬ 
dance matching results use stellar masses assuming a univer¬ 
sal Milky Way IMF (e.g., Ghabrier 2003). A Milky Way IMF 
in spiral galaxies is consistent with dynamical based mass- 
to-light ratios (Courteau & Dutton 2015). However, there 
is increasing evidence for “heavier” IMFs in massive ellipti¬ 
cal galaxies and the bulges of massive spiral galaxies (e.g., 
Gonroy & van Dokkum 2012; Dutton et al. 2013a,b). These 
studies suggest that for massive galaxies (Mstar ^ 10^^M©) 
the observed stellar masses may be underestimated by a 
factor of ~ 2 when assuming a Milky Way IMF and stan¬ 
dard stellar population models. Taking this into considera¬ 
tion nominally improves the agreement of our most massive 
simulations with observations. Note that formally we should 
use the same IMF in the observations and simulations. How¬ 
ever, practically speaking the feedback efficiency parameters 


in the simulations can be simply rescaled. In addition, there 
is enough observational uncertainty in the high mass end of 
the IMF, as well as the actual energy released per SN that 
this inconsistency is not significant. 

The parameters of the feedback model are kept con¬ 
stant across the entire mass range. The parameters were de¬ 
termined based on a comparison with abundance matching 
for one galaxy about the mass of the Milky Way (M 200 ~ 
IO^^Mq) (see Stinson et al. 2013). The agreement between 
the simulated galaxies and the abundance matching at a 
wide range of masses provides strong support for stellar feed¬ 
back being the key ingredient that limits star formation in 
galaxies lower mass than the Milky Way (see also Hopkins 
et al. 2014). 

The points in Fig. [5] are color coded according to the 
number of particles within the virial radius (dark, gas and 
stars). The main haloes have between 4 x 10® and 4 x 10® 
particles. The halo mass range from 2 x 10^® to 10^^ M© 
includes simulations with a wide range of resolutions (from 
20 000 to more than 400 000 particles, and dark matter par¬ 
ticles’ mass is from 6.4 x 10"^ to 1.7 x 10®Mq). There are 
no systematic offsets indicating good convergence in stel¬ 
lar masses at an order of magnitude fewer particles than 
our zoom-in target galaxies. In addition, our central galax¬ 
ies have a range of physical spatial resolutions, with bary- 
onic force softenings varying from 50 pc to 400 pc. This is 
encouraging as cosmological simulations do not always con¬ 
verge when changing the physical resolution (see discussions 
in Vogelsberger et al. 2014; Schaye et al. 2015). 

Below a halo mass of M 200 ~ 2 x IO^^Mq, there starts 
to be a discrepancy between different resolutions. The high 
resolution simulations (blue points) contain less stars than 
lower resolution galaxies. There are two possibilities: Stars 
may form more readily in low resolution galaxies or else 
their halo masses might be reduced through enviromental 
processes such as tidal stripping. While the lower resolution 
halos are outside of the main halos in the z — 0 output, 
it is possible that they have flown nearby the main galaxy 
at some point in the past given their proximity to more 
massive halos. The low resolution halo could form a mass of 
stars according to its earlier halo mass, make a close passage 
of the main halo, and thus have its halo mass reduced while 
retaining its stellar mass (for a study that shows that the 
outermost mass is stripped first, see Chang et al. 2013). 

Above a halo mass of M 200 2 x 10 ^®Mq the scatter 

in the relation is a constant and consistent with observa¬ 
tional constraints of ~ 0.2 dex (More et al. 2011; Reddick 
et al. 2013). Below a halo mass of M 200 ~ 2 x IO'^^Mq the 
scatter starts to increase. More high resolution simulations 
are needed to verify if this feature is due to the relatively 
low resolution of these haloes or environmental effects. 

A number of semi-analytic models have been able 
to produce galaxies that follow the observed relation be¬ 
tween stellar and halo mass at redshift 2 ~ 0 (e.g.. Bower 
et al. 2006; Somerville et al. 2008; Guo et al. 2011). However, 
simultaneously reproducing the evolution has been a chal¬ 
lenge (Weinmann et al. 2012; Hopkins et al. 2014). Fig. [S] 
shows the evolution of the stellar mass vs halo mass rela¬ 
tion since redshift 2 = 4 (a look back time of ~ 12 Gyr). 
Here we only show the most massive halo per simulation. 
Our simulations show good agreement with the abundance 
matching relations from Behroozi et al. (2013) and Moster 
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Figure 5. Stellar mass vs halo mass at redshift 2 = 0 for main simulations (blue points) together with lower mass galaxies in the 
zoom-in region (green and red points). The solid black line and shaded region shows the relation from Kravtsov et al. (2014) derived 
using halo abundance matching. Our simulation matches this very well. The dashed lines show extrapolations of the abundance matching 
relations. For reference, the dotted line shows the cosmic baryon fraction of mass associated with the dark matter halo, indicating that 
our simulations convert only a small fraction of the available gas into stars, as observed. 


et al. (2013). Note that most of the simulated galajdes fall in 
a mass range that can only be extrapolated (dashed lines). 
The range of masses directly constrained by observations is 
shown as the solid lines. 

The agreement with the Mstar — M 200 relation is not 
totally unexpected, as the feedback model employed here 
has been shown to reproduce the time evolution of a sin¬ 
gle Milky Way mass galaxy through the stellar mass vs halo 
mass plane (Stinson et al. 2013), as well as the stellar to halo 
mass relation for galaxies of mass 10®’® < Mstar ^ 10^® ®Mq 
at redshifts z > 2 using lower resolution simulations in a cos¬ 
mological volume of size 100 Mpc (Kannan et al. 2014). 
Our study bridges the gap between these two works by cre¬ 
ating a large number of high-resolution simulations with re¬ 
alistic stellar masses. 

Our simulations are not the first to replicate the Mstar — 
Mhaio relation at 2 = 0 or even at higher redshifts. Aumer 
et al. (2013) presented a sample of 16 galaxies with com¬ 
parable resolution to those presented here but on a smaller 
mass range (10^^ Mq to 2 x 10^^ Mq ). Those simulations 
showed a similar mass and redshift dependence as the sim¬ 


ulations presented in our sample, even though they tend to 
over produce stars at 2 > 4 and under produce them at 
2 < 0.5 in the low mass range, an effect we do not see in our 
(even larger) sample. 

Hopkins et al. (2014) presented a sample of only 7 galax¬ 
ies with a factor of ~ 5 higher mass resolution than ours (see 
Fig- mi, accompanied by several other galaxies that were also 
simulated inside the zoom region. Again, the results are sim¬ 
ilar to what we find, although at all halo masses they over¬ 
produce stellar masses at 2 = 2 when compared with the 
Behroozi et al. (2013) prediction and its extrapolation. All 
three simulation suites were run with different codes and 
stellar feedback recipes, but there are common features to 
all of them. 1) They all include efficient supernova feed¬ 
back (at least 10®^ erg returned to the ISM in all cases). 
2) The supernova feedback is deposited in the local environ¬ 
ment of where the stars formed, promptly, 3 Myr after the 
stars formed. 3) They all include some sort of feedback prior 
to supernova (Aumer: strong, rm ~ 25, radiation pressure; 
Hopkins: radiation pressure, stellar winds, and photoioniza- 
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M200 [Mq] 


Figure 6. Evolution of the Stellar mass vs halo mass relation. Redshift 2 = 0(top left), 2 = l(top right), 2 = 2 (bottom left), and 2 = 4 
(bottom right). Only the most massive halo in each zoom-in region are shown. There is good agreement between our simulations and 
constraints from halo abundance matching (Moster et al. 2013; Behroozi et al. 2013). The dashed lines show extrapolation of the relations 
into mass scales without observational constraints. For reference, the dotted line shows the cosmic baryon fraction of mass associated 
with the dark matter halo. 


tion; this work: strong photoionization included as thermal 
energy). 

Schaye et al. (2015) shows that their simulated stellar 
mass function follows observations at 2 = 0, with an excep¬ 
tion at the knee that translates into a low peak efficiency 
of star formation around Mhaio ~ 10^^ Mq (see their figure 
8). Those simulations invoke a simple prescription for stel¬ 
lar feedback that injects slightly more thermal energy per 
supernova instantaneously after the star particle is formed. 
They add AGN feedback to their energy budget, which pre¬ 
vents over cooling in higher mass galaxies. 

3.2 Star formation rate vs stellar mass 

A consistency check on the stellar mass growth of our simu¬ 
lated galaxies comes from the relation between star forma¬ 
tion rate and stellar mass, sometimes referred to as the “star 
forming main sequence.” 

Fig. [7] shows the evolution of the star forming main 
sequence in our simulations compared with the observed re¬ 
lations (as compiled by Dutton et al. 2010) at 2 ~ 0 (Elbaz 


et al. 2007, Salim et al. 2007), 2 ~ 1 (Elbaz et al. 2007, 
Noeske et al. 2007), 2 ^2 (Daddi et al. 2007) and 2^4 
(Daddi et al. 2009, Stark et al. 2009), together with the 
compilation by Behroozi et al. (2013). Broadly speaking, our 
simulations recover the order of magnitude decline in SFR at 
fixed stellar mass observed since 2 ~ 4.In more detail, they 
match the observations well at all but the lowest mass bins, 
where the observations are more prone to selection biases. 

The agreement at redshift 2 = 2 is especially note¬ 
worthy, as a common feature of galaxy formation models 
is to underpredict the observed SFRs by a factor of > 2 
(e.g., Dutton et al. 2010; Yang et al. 2013; Romeo Velona 
et al. 2013; Somerville & Dave 2014; Furlong et al. 2015). 
This feature of previous models is due to a strong cou¬ 
pling between the cosmological gas accretion and star forma¬ 
tion rates (Weinmann et al. 2011) which leads to SSFR oc 
(1 + In our simulations this coupling is broken by our 

feedback model, which effectively delays star formation at 
early times. Other solutions have included a modihcation to 
the stellar IMF (Dave 2008) or an appeal to systematic bi- 
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Figure 7. Evolution of the galaxy star formation rate vs stellar mass relation. NIHAO simulations are shown with blue points, observed 
relations by lines and shaded regions as indicated. For reference the dashed lines show a specific star formation rate, SSFR=SFR/Matar = 
0.1 Gyr“^. The average SSFR of our simulations decline by a factor of ~ 30 from z = A to the present day. 


ases in the observed star formation rates and stellar masses 
which could plausibly be at the factor of ~ 2 level. 


3.3 Rotation curves 

To this point, our analysis has focused exclusively on the 
mass of stars formed in the simulated haloes. Those com¬ 
parisons show that our simulations are successful according 
to that measure. However, a fully successful galaxy forma¬ 
tion model should also reproduce the mass distribution of 
material in the galaxies. We will look into the matter dis¬ 
tributions in more depth in follow-up papers, but give a 
preliminary look at the rotation curves for the simulated 
galaxies here. 

Fig.[ 8 ]shows the spherically averaged circular velocities 
for all the galaxies in the NIHAO sample. The galaxies have 
been separated into four separate mass groups to limit the 
confusion. The circular velocity curves are colored by their 
halo mass. The curves show a gradual transition based on 
mass that represents a change in the halo response that will 
be explored more explicitly in Dutton et al. {in prep). In gen¬ 


eral terms, the lowest mass halos, M 200 < 2 x IO^^Mq, have 
rotation curves that rise relatively steeply. As the masses 
increase, the rotation curves rise more slowly. Then, the 
galaxies cross a higher mass threshold around 6 x lO^^M© at 
which the galaxies revert to a steeply rising rotation curve, 
which in some cases grows to have a high central peak. Such 
centrally peaked rotation curves do exist in observed mas¬ 
sive spiral galaxies (e.g., Noordermeer et al. 2007; Dutton 
et al. 2013a). However, we note that our simulations do not 
include any form of AGN feedback. The central peaks in the 
rotation curves may thus represent over cooling that would 
otherwise be heated by AGN or another higher energy phe¬ 
nomenon. 


4 SUMMARY 

We introduce a new, unmatched sample of high resolu¬ 
tion hydrodynamic simulations from the NIHAO project. 
The galaxies are individually simulated in their own zoom 
regions. They span a wide range of halo masses, 5 x 
10® <M 2 oo/M 0 < 3 X 10*^® and stellar masses 5 x 10"* < 
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Figure 8 . Rotation curves of the spherically averaged circular velocity GM/r, split into four halo mass ranges, as indicated. 


Mstar/M 0 < 2 X 10^^. Across this range of masses, the galax¬ 
ies show a range of morphologies, colors and sizes that corre¬ 
spond well with observed galaxies. The only energetic source 
term for the simulations comes from stars; stellar feedback 
due to photoionization from hot, young stars; and the ther¬ 
mal energy from supernovae. The stellar masses of the sim¬ 
ulations closely track the results from abundance matching 
for stellar masses as a function of halo mass across more 
than three quarters of the life of the Universe. That means 
that a complete cosmological sample of these galaxies would 
populate the stellar mass function as observed. The corre¬ 
spondence of the simulations to observations provides a final 
proof that stellar feedback is the chief piece of physics re¬ 
quired to limit the efficiency of star formation in galaxies 
less massive than the Milky Way. 
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APPENDIX A: PARAMETERS 

In this appendix we summarize all the galaxy and halo pa¬ 
rameters for main target haloes at redshift z = 0 through 
this paper. Table |XT] contains the simulation ID, the number 
of total particles, dark matter particles and star particles, 
virial mass, stellar mass and star formation rate. 
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Table Al. Galaxy and halo parameters for main target haloes at redshift 2 = 0 


simulation ID 

^ particles 

?^200 

# DM particles 

>^DM 

# star particles 

^star 

virial mass 
(M' 2 Oo/M 0 ) 

stellar mass 

(A^star/M0) 

star formation rate 
(M 0 yr-1) 

g3.54e09 

1,200,374 

1,157,470 

85 

3.95 

X 

10 ® 

3.02 

X 

lO'^ 

0.00 



g4.36e09 

1,278,589 

1,167,760 

79 

9.54 

X 

10 ® 

3.30 

X 

lO'^ 

0.00 



g4.99e09 

320,446 

295,765 

538 

5.72 

X 

10 ® 

3.76 

X 

10 ® 

5.8 X 

10 “ 

-5 

g5.22e09 

353,754 

327,288 

173 

6.32 

X 

10 ® 

1.18 

X 

10 ® 

1.2 X 

10 “ 

-5 

g5.41e09 

466,210 

443,516 

2,956 

5.11 

X 

10 ® 

1.21 

X 

10 ® 

0.00 



g5.59e09 

358,332 

335,527 

2,520 

6.46 

X 

10 ® 

1.74 

X 

10 ® 

1.7 X 

10 “ 

-4 

g5.84e09 

555,955 

512,178 

87 

5.95 

X 

10 ® 

1.16 

X 

10 ® 

6.3 X 

10 “ 

-5 

g7.05e09 

584,878 

544,306 

3,249 

1.05 

X 

lOlo 

2.27 

X 

10 ® 

0.00 



g7.34e09 

559,724 

521,731 

331 

6.05 

X 

10 ® 

4.22 

X 

10 ® 

0.00 



g9.26e09 

581,361 

527,275 

132 

6.14 

X 

10® 

5.37 

X 

lO'^ 

0.00 



gl.09el0 

690,285 

546,318 

9,204 

1.09 

X 

lOlo 

6.55 

X 

10® 

0.00 



gl.lSelO 

628,317 

565,451 

4,944 

1.10 

X 

lOlo 

3.41 

X 

10® 

0.00 



gl.23el0 

515,485 

468,423 

2,314 

9.08 

X 

10® 

1.60 

X 

10® 

0.00 



gl.44el0 

1,159,093 

834,135 

9,656 

1.70 

X 

lOio 

6.69 

X 

10® 

0.00 



gl.47el0 

862,592 

784,867 

13,043 

1.52 

X 

lOio 

9.00 

X 

10® 

1.3 X 

10- 

-3 

gl.50el0 

791,452 

725,701 

5,055 

1.40 

X 

lOio 

3.42 

X 

10® 

0.00 



gl.57el0 

728,314 

670,688 

13,027 

1.29 

X 

lOio 

8.93 

X 

10® 

0.00 



gl. 88 el 0 

1,063,823 

971,852 

23,864 

1.88 

X 

lOio 

1.65 

X 

10" 

3.9 X 

10“ 

-3 

gl.89el0 

1,401,592 

1,060,626 

17,813 

2.13 

X 

lOio 

1.29 

X 

10" 

3.1 X 

10- 

-3 

gl.90el0 

1,457,806 

1,222,038 

17,279 

2.40 

X 

lOlo 

1.20 

X 

10" 

1.4 X 

10- 

-3 

gl.92el0 

1,272,888 

990,462 

8,255 

1.98 

X 

lOlo 

5.30 

X 

10® 

1.4 X 

10- 

-4 

gl.95el0 

794,621 

704,686 

5,488 

1.37 

X 

lOlo 

3.82 

X 

10® 

0.00 



g2.09el0 

1,018,898 

839,068 

13,600 

1.66 

X 

lOlo 

9.35 

X 

10® 

0.00 



g2.34el0 

1,500,299 

1,319,369 

19,729 

2.57 

X 

lOlo 

1.36 

X 

10" 

0.00 



g2.37el0 

1,503,953 

1,232,912 

25,627 

2.43 

X 

lOlo 

9.65 

X 

10® 

0.00 



g2.39el0 

869,294 

743,922 

8,601 

1.46 

X 

lOlo 

5.94 

X 

10® 

1.7 X 

10- 

-4 

g2.63el0 

458,723 

414,291 

18,388 

2.70 

X 

lOio 

4.28 

X 

10" 

0.00 



g2.64el0 

2,340,543 

1,579,110 

43,275 

3.26 

X 

lOio 

2.93 

X 

10" 

0.00 



g2.80el0 

521,957 

413,400 

15,628 

2.77 

X 

lOio 

3.68 

X 

10" 

2.5 X 

10- 

-3 

g2.83el0 

483,150 

370,911 

13,428 

2.50 

X 

lOio 

2.96 

X 

10" 

0.00 



g2.94el0 

1,908,993 

1,657,184 

85,637 

3.22 

X 

lOio 

5.86 

X 

10" 

4.2 X 

10- 

-4 

g3.19el0 

677,287 

525,030 

10,679 

3.54 

X 

lOio 

1.47 

X 

10" 

3.9 X 

10- 

-5 

g3.44el0 

1,077,981 

694,131 

26,561 

4.88 

X 

lOio 

6.32 

X 

10" 

4.7 X 

10- 

-3 

g3.67el0 

547,857 

480,668 

23,584 

3.15 

X 

lOio 

5.49 

X 

10" 

0.00 



g3.93el0 

682,659 

477,948 

16,007 

3.30 

X 

lOio 

3.75 

X 

10" 

0.00 



g4.27el0 

885,293 

616,773 

26,191 

4.25 

X 

lOio 

6.16 

X 

10" 

4.4 X 

10- 

-3 

g4.48el0 

1,191,821 

895,736 

58,334 

6.04 

X 

lOio 

1.37 

X 

10* 

3.9 X 

10- 

-5 

g4.86el0 

1,051,559 

757,795 

52,672 

5.16 

X 

lOio 

1.22 

X 

10* 

3.7 X 

10- 

-3 

g4.94el0 

984,362 

791,322 

48,035 

5.26 

X 

lOio 

1.11 

X 

10* 

7.8 X 

10- 

-3 

g4.99el0 

940,085 

732,808 

52,511 

4.90 

X 

lOio 

1.22 

X 

10* 

9.4 X 

10- 

-3 

gS.OSelO 

777,412 

650,173 

40,840 

4.29 

X 

lOio 

9.47 

X 

10" 

7.8 X 

10- 

-5 

g 6 . 12 el 0 

986,798 

733,881 

39,024 

4.97 

X 

lOio 

9.13 

X 

10" 

7.8 X 

10- 

-5 

g6.37el0 

2,621,929 

1,617,864 

95,020 

1.15 

X 

loll 

2.12 

X 

10* 

1.2 X 

10- 

-2 

g6.77el0 

2,118,102 

1,333,784 

204,017 

9.28 

X 

lOio 

4.84 

X 

10* 

8.7 X 

10- 

-2 

g6.91el0 

1,331,679 

1,069,321 

107,120 

7.08 

X 

lOio 

2.50 

X 

10* 

2.9 X 

10- 

-3 

g7.12el0 

1,009,755 

715,315 

58,766 

4.88 

X 

lOio 

1.37 

X 

10* 

1.5 X 

10- 

-2 

g8.89el0 

591,924 

397,347 

50,174 

9.22 

X 

lOio 

4.02 

X 

10* 

2.9 X 

10- 

-2 

g9.59el0 

614,423 

368,306 

34,468 

8.84 

X 

lOio 

2.76 

X 

10* 

8.7 X 

10- 

-2 

gl.OSell 

775,657 

505,507 

70,478 

1.18 

X 

loll 

5.67 

X 

10* 

6.0 X 

10- 

-2 

gl.OSell 

785,537 

520,108 

108,025 

1.20 

X 

loll 

8.47 

X 

10* 

2.7 X 

10- 

-2 

gl.37ell 

1,014,522 

653,073 

252,530 

1.48 

X 

loll 

2.02 

X 

10® 

9.5 X 

10- 

-2 

gl.52ell 

1,128,519 

657,385 

112,285 

1.57 

X 

loll 

7.92 

X 

10* 

3.8 X 

10- 

-2 

gl.57ell 

1,076,183 

677,124 

143,921 

1.58 

X 

loll 

1.15 

X 

10® 

6.7 X 

10- 

-2 

gl.59ell 

1,223,754 

691,453 

86,944 

1.68 

X 

loll 

6.69 

X 

10* 

7.9 X 

10- 

-4 

gl.64ell 

1,468,563 

786,245 

111,008 

1.93 

X 

loll 

9.13 

X 

10* 

0.38 



g2.04ell 

1,722,119 

901,588 

573,355 

2.09 

X 

loll 

4.70 

X 

10® 

0.79 



g2.19ell 

920,447 

557,247 

113,958 

1.31 

X 

loll 

9.27 

X 

10* 

9.8 X 

10- 

-2 

g2.39ell 

2,133,166 

1,109,605 

698,696 

2.58 

X 

loll 

5.80 

X 

10® 

1.12 



g2.41ell 

1,939,348 

1,092,825 

501,991 

2.53 

X 

loll 

4.10 

X 

10® 

0.62 



g2.42ell 

2,073,620 

1,175,947 

689,091 

2.68 

X 

loll 

5.48 

X 

10® 

0.39 



g2.54ell 

2,004,595 

1,144,673 

424,390 

2.67 

X 

loll 

3.50 

X 

10® 

0.87 
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Table A2. cont’d Galaxy and halo parameters for main target haloes at redshift z = 0 


simulation ID 

^ particles DM particles 

^^200 

# star particles 

^star 

virial mass 
(M'2 Oo/M0 ) 

stellar mass 

(A^star/M0) 

star formation rate 
(M0 yr-i) 

g3.06ell 

2,688,990 

1,323,580 

906,552 

3.11 

X 

IQll 

7.42 

X 

10 ® 

1.25 

g3.21ell 

2,342,924 

1,272,210 

448,644 

3.03 

X 

IQll 

3.67 

X 

10 ® 

1.10 

g3.23ell 

631,510 

371,614 

44,745 

8.93 

X 

IQlO 

3.60 

X 

10 * 

5.5 X 10"® 

g3.49ell 

3,276,200 

1,815,975 

539,486 

4.33 

X 

IQll 

3.97 

X 

10 ® 

0.61 

g3.55ell 

3,241,957 

1,760,597 

480,380 

4.23 

X 

IQll 

3.85 

X 

10 ® 

0.66 

g3.59ell 

2,762,536 

1,459,228 

534,349 

3.49 

X 

IQll 

4.36 

X 

10 ® 

0.74 

g3.71ell 

4,015,748 

1,691,651 

1,535,547 

4.08 

X 

IQll 

1.24 

X 

lOio 

1.49 

g4.90ell 

2,452,320 

1,369,494 

428,202 

3.25 

X 

IQll 

3.43 

X 

10 ® 

0.42 

g5.02ell 

5,258,849 

2,408,247 

1,821,403 

5.75 

X 

IQll 

1.46 

X 

lOio 

0.92 

gS.Slell 

3,534,992 

2,204,411 

559,381 

5.28 

X 

IQll 

1.64 

X 

lOio 

2.26 

g5.36ell 

4,505,028 

2,889,033 

445,661 

6.88 

X 

IQll 

1.18 

X 

lOlo 

3.88 

g5.38ell 

4,187,358 

2,726,961 

648,971 

6.46 

X 

IQll 

1.85 

X 

lOio 

1.24 

g5.46ell 

2,459,236 

1,374,445 

470,795 

3.25 

X 

IQll 

3.78 

X 

10 ® 

0.37 

g5.55ell 

5,020,346 

2,200,598 

2,075,884 

5.20 

X 

IQll 

1.71 

X 

lOio 

2.75 

g6.96ell 

1,094,722 

412,962 

498,881 

8.00 

X 

IQll 

3.30 

X 

lOio 

10.40 

g7.08ell 

1,098,831 

412,287 

479,668 

8.06 

X 

IQll 

3.05 

X 

lOio 

3.66 

g7.44ell 

1,294,128 

596,875 

283,637 

1.18 

X 

IOI2 

1.85 

X 

lOio 

12.95 

g7.55ell 

1,185,149 

455,930 

483,752 

8.93 

X 

IQll 

3.11 

X 

lOio 

2.92 

g7.66ell 

1,528,867 

478,007 

911,988 

9.30 

X 

IQll 

5.92 

X 

lOio 

2.71 

g8.06ell 

1,366,038 

481,349 

665,796 

9.43 

X 

IQll 

4.48 

X 

lOio 

13.19 

g8.13ell 

1,704,130 

505,171 

1,043,618 

9.91 

X 

IQll 

6.68 

X 

lOio 

3.67 

g8.26ell 

1,513,265 

518,236 

739,749 

1.02 

X 

IOI2 

4.68 

X 

lOio 

2.04 

g8.28ell 

1,276,937 

591,330 

275,702 

1.16 

X 

IOI2 

1.77 

X 

lOio 

1.01 

gl.l2el2 

1,977,112 

564,785 

1,222,359 

1.12 

X 

IOI2 

7.85 

X 

lOio 

4.29 

gl.77el2 

3,668,318 

1,117,439 

2,173,715 

2.19 

X 

IOI2 

1.39 

X 

loll 

8.47 

gl.92el2 

4,018,048 

1,200,667 

2,467,821 

2.34 

X 

IOI2 

1.58 

X 

loll 

8.13 

g2.79el2 

5,598,386 

1,800,095 

3,099,962 

3.53 

X 

IOI2 

1.96 

X 

loll 

11.09 
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